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Abstract 

QCD with eight flavors is studied on 16^ x Nf lattices with Nt = 4, 
Oh! 6, 8, 16 and 32, a dynamical quark mass ma = 0.015 and lattice cou- 

pling (3 = 6/g'^ between 4.5 and 5.0. For A^^ = 16 and 32, hadron 
masses and screening lengths are computed for a variety of valence 
quark masses. The previously observed, strong, first-order transition 
kS ' for A^t = 4, 6 and 8 is seen, for Nt = 16, to become a /3- independent, 

C^ ' zero-temperature transition characterized by a factor of ~ 3 change 

in lattice scale. This strong, first-order transition restores chiral sym- 
metry, at least for A^f = 4, 6 and 8, producing a chirally symmetric, 
weak-coupling phase. However, as Nt increases to 16, the chiral sym- 
metry properties of the weak-coupling side of the zero-temperature 
transition are unclear and offer a hint of a normal, finite-temperature, 
chiral symmetry breaking transition in the weak-coupling phase. 

PACS numbers: 12.38.Gc, 11.15.Ha 
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1 Introduction 

Among the possibilities offered by the numerical simulation of Quantum 
Chromodynamics (QCD) is that of studying a variety of values of the phys- 
ical parameters of the system. In particular, a deeper understanding of the 
physics of QCD may result from studying its dependence on the dimension 
of the fundamental representation of the gauge group (the number of colors, 
Nc) and the number of quark flavors, Nf. In this paper we pursue the later 
alternative by considering the case of QCD with eight light quark flavors. 

By considering a range of lattice shapes and couplings we intended to 
study the effects of eight quark flavors on both the finite-temperature QCD 
phase transition as well as the hadron spectrum at zero temperature. How- 
ever, as we will see, such an investigation of the hadron spectrum is seriously 
impeded by the complex phase structure of eight-flavor QCD. 

Our work extends earlier T > 0, eight-flavor calculations [||, ||, ||, ^ to 
larger lattices, smaller quark masses and smaller lattice spacings, with per- 
haps surprising results. The earlier work on 8^ x 4[jl|, 0], 6^0, 8^0, §] and 
16^x4 and x6[Q lattices shows a quite strong, first-order transition. Compar- 
ing the transition for Nf = 2, 3, 4 and these eight-flavor calculations reveals a 
strengthening of the transition as the number of flavors increases|^, |^. Both 
greater metastability at the critical value of (3 and an expanded window of 
small quark mass within which the transition occurs are found as Nf is in- 
creased. Regardless of the number of flavors, significant variation of (3c is 
seen as the number of sites in the temperature direction, Nt, is increased 
from 4 to 6, as is expected for a "finite-temperature" transition. 

As we show below, this pattern changes significantly for larger lattices 
with Nf = 8. First, the variation of (3c with Nf has vanished for A^^ > 8; 
both 16^ X 8 and 16^ x 16 lattices show a strong, first-order transition at 
the same value of (3. Thus the temperature dependent transition seen for 
Nt < 8 has become a temperature independent "bulk" transition for Nf > 8, 
which we describe as separating "strong"- and "weak"- coupling phases of 
the eight-flavor theory. (The presence of a bulk transition for large A^t is 
suggested by earlier At-independent jumps seen in Nf = 8 calculations with 
heavier quarks. 0) Second, although the valence quark mass dependence of 
(xx) suggests that the transition seen for A^ = 4, 6 and 8 is one of chiral- 
symmetry restoration (as is the case for Nf = 2, 3 and 4), the situation is less 
clear for Nt = 16, where (xx) shows non-linear behavior for small valence 



quark mass in the weak-coupling phase. 

However, if /3 is increased to ~ 0.3 above f3c for the bulk transition, chirally 
symmetric behavior is seen for both (xx); which now depends linearly on the 
valence quark mass, and the hadron correlation functions. This suggests that 
in addition to the first-order, bulk transition separating strong- and weak- 
coupling phases, the weak-coupling phase may itself be divided into two 
distinct phases separated by a normal, finite-temperature (i.e., Nt dependent) 
transition or cross-over region. In fact, the non- linear behavior seen for (xx) 
in the weak-coupling phase is reminiscent of the mass dependence of (xx) 
seen in four-flavor calculations on 10^ x 6 when /3 ?^ /5c0. 

This behavior is described by the phase structure in the j3 — Nt plane 
shown in Figure |I]. (The figure depicts a lattice of infinite spatial extent.) 
The solid line represents the strong, first-order transition, which varies with 
Nt for Nt < 8 and becomes a zero-temperature or bulk transition for A^^ > 8. 
The dashed line expresses our speculation that a normal, finite-temperature 
transition or cross-over region is also present in the weak-coupling phase. 
This line is drawn with a slope given by the perturbative renormalization 
group, since this gives the dependence of /3c on Nt as Nt —^ oo for a physical 
transition/cross-over. The dotted line is an extension of the dashed line to 
values of (3 smaller than /3c for the bulk transition, since the strength of 
the bulk transition and the properties of our updating algorithm make it 
possible to study the weak-coupling phase when it is only metastable. The 
solid squares locate parameter values that we have studied and the open 
squares mark values of /3c that have been identified. 

Unfortunately, a study of T = hadron masses in the weak-coupling 
region must be performed within the wedge-shaped region in Figure |l| that 
is bounded below by the dashed line and to the left by the solid line. As 
will be discussed, volumes larger than 16^ will be required to unambiguously 
recognize this region, let alone perform a meaningful mass calculation there. 

In Section |^ we describe the calculations that have been performed and 
the methods used in both generating the gauge configurations studied and 
constructing the various observables computed. The transition for iVt = 4, 6 
and 8 is considered in Section H, where we present evidence that it separates 
chirally symmetric and asymmetric phases. In Section | the bulk transition, 
isolated on 16^ x 16 and 16'^ x 32 lattices, is discussed while in Section || we 
consider the weak-coupling phase on 16^ x 32 volumes for larger (3 and hence 
higher temperature. In Section ^ we discuss the relation between the phase 



structure presented here for eight flavors with that seen for smaller numbers of 
flavors. We suggest that the familiar cross-over region separating strong and 
weak coupling in the zero-flavor theory strengthens with increasing Nf and 
becomes our Nf = 8 bulk transition. Finally in Section |^ various concluding 
remarks and speculations are presented. 



2 Description of the Calculation 

We have carried out a Monte Carlo evaluation of the Euclidean-space, Feyn- 
man path integral for full QCD using the R algorithm of Gottlieb, et al.^. 
The calculation required about five months on the 256-node Columbia ma- 
chine, a 16 X 16 mesh of fast array-processors which achieves a sustained 
performance for these calculations of 6.4 Gflops^]. The R algorithm evolves 
the gauge fields according to the action 

S = -^f3Y.RetTUr- jNf In det[{D + ma) [D^ + ma)]. (1) 

Here the first term is the usual Wilson action with U-p the product of the 
four SU{3) link matrices that border the plaquette V. The second term rep- 
resents the effects of Nf degenerate flavors of dynamical fermions of mass m. 
The factor of 1/4 preceding this term compensates for the fermion doubling 
present in the staggered Dirac operator D. The additional doubling intro- 
duced in Eq. (0) by squaring the Dirac operator is removed by restricting the 
squared operator to even lattice sites. The operator D can be defined by its 
action on an S'f/(3)-triplet field 0: 

Here Un,^i and 77„,^ are the link matrix and staggered fermion sign factor 
associated with a lattice link extending from the site n in the direction /i to 
the site n + /x. 

In simulating eight flavors of staggered fermions, we have a choice be- 
tween two well established methods: the R algorithm, which contains finite 
time-step errors of order (Ar)^, and the exact Hybrid Monte Carlo method 
of Duane, et a/.f^, which requires twice the number of Dirac propagator in- 
versions per unit of Monte Carlo time when used to simulate Nf = 8. We 



chose the R algorithm with a time step Ar = 0.0078125 for this exploratory 
calculation both because code to perform the exact, eight-flavor update was 
not ready when we wished to begin the simulation and we wanted to re- 
duce the required computer time. We have explicitly studied the effects of 
the finite time-step errors and made a comparison with results from the Hy- 
brid Monte Carlo algorithm in our determination of (3c for the transition on 
16^ X 8 and 16^ lattices as is discussed in Sections | and ^. The comparisons 
presented there show the expected quadratic dependence on Ar. Although 
the Ar-errors found are quantitatively large {e.g. 5% in (3c) , the qualitative 
features of the calculation appear unaffected. 

Thus, except where otherwise noted, this calculation is performed with 
a time step Ar = 0.0078125 and a molecular dynamics trajectory of length 
0.5 time units. We have used three types of starting configurations in this 
calculation: hot starts where the gauge fields are disordered, cold starts 
where all the gauge link matrices are unit matrices and mixed starts which 
are described in detail in Section ^. After each trajectory, the molecular 
dynamics "momenta" are randomized and measurements on the link variables 
carried out. In particular, after each trajectory we compute average values 
of the Wilson action and the fermion operator {xx)- Our gauge action is 
(1 — 1/3-Retr f/p) and our convention for (xx) is 

{xx) = 3 j;pj;^ Y. (^«^") ' (3) 

where the sum is over all points in the lattice, (xx) is estimated by 

(xx) = lNh{{^l,nhl{-l^)^^^hn)) 

where for each site n, /i„ is an independent, complex three- vector of Gaussian 
random numbers and ((■ ■ ■)) denotes an average over gauge fields and the 
random three- vectors, /i„. For the work in this paper, we have used three 
sets of hnS for each gauge configuration. Also, we restrict the h^s to even 
sites in evaluating the squared operator in Eq. (|) and multiply the result by 
two. 

Hadron propagators are calculated every 5 units of microcanonical time 
from quark propagators determined using Coulomb gauge wall sources whose 



spatial size is the spatial volume of the lattice. The quark propagator from 
a source at time slice t is calculated, for each color index a = 1-3, using as a 
source an SU{3) triplet field /i^j/, given by 

where n = n mod 2. For a given source time slice, the three quark prop- 
agators are then combined into hadron propagators corresponding to ma- 
trix elements of the five conventional local hadron operators |]lOl- In order 
to improve the statistical accuracy of our results, we use the average-over- 



time-slice (AOTS) method|Tl|l in which this hadron propagator calculation is 
performed A^^ times placing the wall source on each time slice in the lattice. 
The resulting Nt propagators are then averaged together. 

In both the updating steps and the calculation of the hadron propagators 
we must solve a Dirac equation of the form {D + ma)y = h. We perform the 
required inversion oi D + ma using the conjugate gradient algorithm. We 
iterate this method until our approximate solution after the ith iteration, yi, 
yields an appropriately small residual vector r^ = {D^D + {maY)yi — {D + 
ma)h. Specifically, we perform the inversion on the even sub lattice (since 
the solution on the odd sublattice can be found from this) and iterate until 



^{r„r,)/{h,h) < A, (6) 

where the inner product (a, b) of the complex vectors is over even lattice 
sites and colors. For the inversions that occur in the updating steps and 
{xx) "we use A = 6.38 x 10~^/v^iVt and perform typically between 300 and 
700 conjugate gradient iterations, depending on (3, for At = 32. For the 
hadron propagator calculation we use the somewhat more stringent condition 
A = 2.21 X 10-^ yielding 700 to 800 iterations for Nt = 32 with ma = 0.015. 
Much is to be learned in calculations of this sort by varying the quark 
mass used in the simulation. Perhaps of greatest interest is the variation of 
(xx) ^-nd m^ with quark mass. A non-zero value of (xx) ^^nd a zero value 
ttItt in the m — *• limit are both definitive indicators of the spontaneous 
breaking of chiral symmetry. Unfortunately, in the exploratory calculation 
reported here only the single value m,a = 0.015 has been used. However, 
we have computed the dependence of {xx) ^-nd the hadron masses on the 
quark mass that appears in the quark propagators that explicitly enter the 
evaluation of the right-hand side of Eq. (^ and the hadron masses. We have 
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only considered the case where all quark propagators used to form a hadron 
propagator have the same quark mass. 

Thus we distinguish two quark masses that enter our calculation: the 
"sea" quark mass rusea that enters the quark determinant in the path inte- 
gral and the "valence" mass rriyai which appears in the quark propagators 
that make up the various observables. In this way we can define, for exam- 
ple, m.„{msea-i^nn^ai)- A consistent calculation with Nf fiavors of degenerate 
quarks requires rrivai = in sea- We might call the quantity m^(msea,m„a/) a 
"quenched" approximation to the proper quantity mT.imyaU'i^vai)- However, 
a non-vanishing limit of (xx) ^is iriyai ^ is never-the-less an indicator of 
spontaneous symmetry breakdown, although the observable in question is no 
longer local. Likewise the Goldstone theorem implies that if this quenched 
{'X.x)ij^seai i^vai) IS uou- Vanishing in the limit m^ai -^ 0, then the correspond- 
ing mTj{msea,mvai) will also vanish in that limit. Clearly, the limit rriyai — >■ 
provides us with interesting information about the character of the small 
eigenvalues of the Dirac operator. 

In fact, if the masses rrisea and rrivai are sufficiently small that m^ and 
(xx) depend on them linearly, and if (xx) is non-vanishing for both the limits 
f^vai = T^sea — > and rUsea fixcd, m^ai — ^ 0, then necessarily the two values 
of m^ agree: 

ml {rrisea , m^ai ) = ml {m^ai , m^ai ) +0{m^), (7) 

since both sides are linear functions of m^ai which agree at two points, m^ai = 
rUsea and rriyai = 0. Because the majority of the linear mass dependence of 
{xx) comes from the quadratically divergent term proportional to rriyai we 
might also expect this quenched calculation of {xx) to be quite accurate. Our 
Nf = 2 calculations [|I^], in which four values of nisea were used, bear these 
expectations out. In Figure |] we show ml and {xx) fo^' the four normal 



points|]l2| m^aiO' = rngeaO- = 0.01, 0.015, 0.02, 0.025 and a fifth, quenched 
point||l^ m^seaCi = 0.01, m.^aiO' = 0.004. The values for the quenched point 
are m^^^a = 0.173(5) and {xx) = 0.0157(1). A linear fit to these five points 
yields: 



2^2 



mia 



0.0033(15) + 5.86(7) m„,da {x^ dof = 26/3), 



{xx)a^ =0.00768(11) + 2.034(6) m„,rfa {x^/dof = U/3). ^^^ 

Although the x^ ior these fits is large, one can see from Figure H that the 
m^sead = 0.015 point is the dominant contributor to the large x^, not the point 
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at ma = 0.004. Leaving out the rriseaO' = 0.015 point gives fits with slopes and 
intercepts that are the same within errors and which have x^ /dof = 0.2/2 
and 0.5/2, respectively. 

3 Chiral-symmetry Restoration: A/^ = 4, 6 
and 8 

As noted in the introduction, previous studies on lattices with Nt < 8 have 
observed a strengthening of the chiral transition for an increasing number of 
quark flavors. For Nf = S, (xx) has been seen to change by about a factor 
of two across the transition for quark masses of 0.1 on lattices with Nt = 4, 
6 and 8|0, Q. Previous eight-flavor work did not include an extrapolation 
of (xx) to zero quark mass in the weak-coupling phase to demonstrate that 
the transition restored chiral symmetry. In this section we report on our 
investigation of this transition for Nt = A, 6 and 8 lattices — both the accurate 
determination of Pc and the m^ai —^ extrapolation of {xx) which gives 
evidence that this transition does restore chiral symmetry. Our A''^ = 4 
results were obtained from a series of simulations described in Table | while 
Tables |T| and |ITT| contain a similar description of the Nt = 6 and 8 runs. All 
the lattices had a spatial volume of 16^. 

The accurate determination of /3c for a strong, first-order transition presents 
a familiar dilemma: if we work with a large spatial volume, the considerable 
metastability of both phases implies a large range of (3 within which each 
phase appears to be stable, even for quite long Monte Carlo evolution times. 
(At least if we use existing, local updating algorithms.) If a sufficiently small 
spatial volume is used to eliminate this metastability, the transition may 
be significantly distorted by finite-volume effects. We solve this problem by 
working with a large volume but beginning with a configuration in a mixture 
of phases. Starting with such a mixed phase, very small changes in the choice 
of f3 cause the system to rapidly evolve into either of the two phases |T^ . 



Our results for f3c and {xx) at At = 4, 6 and 8 have been obtained by 
starting from mixed-phase configurations generated as follows: An initial /3 
value was chosen for which hot and cold starts gave two metastable phases, 
with different values of (xx) and gauge action, which were stable for 50 
or more units of microcanonical time. We then chose a configuration from 



one phase and evolved it, changing (3 every 10-20 units of time until the 
configuration had a value of (xx) ^"^^ action close to halfway between their 
values in the two phases. (Both (xx) ^ind the action reached their halfway 
values concurrently.) The final value of j3 used in evolving this mixed-phase 
configuration was chosen to keep both (xx) ^ind the action roughly constant 
for 10-20 time units, thus diminishing any "inertia" the lattice might have 
pushing it in the direction of either phase. 

The upper curve in Figure shows the evolution of (xx) produced when 
generating a mixed-phase configuration from a hot start for an A^t = 4 lattice. 
The first 80 time units show normal evolution with (3 = 4.5. During the next 
95 time units /? was continually adjusted to produce the intermediate value 
of (xx) shown. 

After generating a mixed-phase configuration, /?c was found by perform- 
ing a series of evolutions, with different values of /3, each starting from the 
given mixed-phase configuration. An example of this procedure is repre- 
sented in Figure §. There we show the evolution of (xx) for five values of 
/3, starting from the mixed-phase configuration generated from an initial hot 
start. Clearly (3 = 4.5 lies on the strong-coupling side of the transition while 
j3 = 4.65 falls on the weak-coupling side. The slow evolution of the /? = 4.58 
run locates Pc while the short runs at 4.55 and 4.6 suggest Pc = 4.58(1) as a 
reasonable conclusion for f3c with errors. 

In order to demonstrate the reliability of this procedure we repeated the 
determination using a second, mixed configuration generated instead from 
an initial cold start. The evolutions for three choices of f3 shown in Figure ^ 
behave in a manner very consistent with Figure]^ and the value /3c = 4.58(1) 
deduced above. This indicates that our mixed-phase configurations are in- 
dependent of whether they were made from a hot or cold start and argues 
against any bias toward one phase. In addition, we have never seen an evo- 
lution beginning from one of our mixed-phase configurations begin to change 
in the direction of one phase and then reverse itself. This indicates that our 
mixed-phase configurations have no inertia toward a particular phase. 

The series of evolutions used to determine (3c for A^f = 8 is shown in 
Figure H, from which we deduce /?c = 4.73(1). Although not shown, the 



evolutions of the action are very similar to those for (xx) • Table IV lists our 
results for (3^ for A^t = 4, 6 and 8. For comparison this table also includes 
the Nt = 16 results which are discussed in the next section. 

Given the ease with which we can determine (3c and its importance in our 
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later considerations, it is reasonable to study the effects of finite time-step 
errors by computing f3c for a number of choices for Ar. We concentrated on 
the Nt = 8 case with the results given in Table |rV|. As is indicated in Table |T| 
these values of Pc were determined by the procedure described above, starting 
from the same mixed phase (generated with Ar = 0.0078125) and evolving 
using a series of updating schemes with Ar = 0.002, 0.005, 0.0078125 and 
0.0125 and with the exact Hybrid Monte Carlo algorithm. 

Figure ^ shows the dependence of Pc on Ar. Leaving out the point at 
Ar = 0.0125 and the point from the exact algorithm, we find that (3c is fit 
by 

/5c(Ar) = 4.58(1) + 2460(250) (Ar)^ (9) 

with x^/dof = 0.02/1. (The fit for (3^ using the points up to 0.0078125 is 
good enough to make it hard to fit the 0.0125 point, even if higher order 
terms are included.) Clearly, we find the expected (Ar)^ dependence of Pc 
and good agreement between the constant term in the fit and the value of Pc 
from the exact algorithm. 

Finally, let us examine the values of (xx) obtained from these various 
runs and their dependence on m^ai- Tables 0, ^ and |VII| give (xx) ^or 



the masses used. The thermalization times, rgg, given in these three tables 
are estimates obtained by eye from the plots of (xx) and vary because the 
thermalization time depends on P — P^, which is not constant for the different 
runs. 

As can be seen in Figure ^, {xx) extrapolates linearly to zero as m^ai -^ 
for iVt = 4, 6 and 8 on the P > Pc side of the transition. The fits are forced 
through the origin and have x'^/dof = 0.39/1, 3.9/2 and 5.1/1, respectively. 
The figure clearly shows that our results are consistent with a chirally sym- 
metric, weak-coupling phase. In addition, for A^f = 8 we can check that 
similar, chirally symmetric behavior is seen as the time step is varied by ex- 
amining Ar = 0.0125 and 0.005. For these cases we also find good linear 
fits for (xx) a-s a function of quark mass. When forced through zero, the fits 
have x'^/dof = 0.74/1 and 0.007/1, respectively. 

4 Evidence for a T=0 transition 

The eight-flavor results described in the preceding section and earlier work 
of others look much like the chiral-symmetry-restoring phase transition seen 
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for four flavors of light quarks on lattices of similar size. These results are 
usually interpreted as a lattice approximation to a phase transition at non- 
zero temperature for the continuum field theory. The variation of the value 
of (3 where the transition occurs (/3c) with Nt supports this interpretation. 
One expects the physical temperature of the lattice at the critical point, 
T = {Nta)~^, to be fixed so that changes in the lattice spacing a resulting 
from changes in (3 must be compensated by changes in A^^. Such behavior 
is quite well established in pure QCD where the variation of a (and hence 



Nt) with P predicted by the perturbative renormalization group is seen[p!5 
on lattices as large as 24^ x 16. 

However, for the eight-flavor transition we find the critical value of (3 
does not change when Nt is increased from 8 to 16. For Nt = 16 we continue 
to see a very strong, first-order transition even though the 16^ lattice now 
has a spatial size no larger than the temporal extent. This apparent space- 
time volume independence of Pc suggests the transition will persist with this 
fixed value of jSc even for a system of infinite spatial and temporal extent. We 
conclude that for Nf = 8 there is a strong, first-order T = or bulk transition 
separating the strong- (/3 < Pc) and weak- (/? > Pc) coupling regimes. 

Although the strong-coupling phase seen on 16^ lattices appears much 
like that found for iVt = 4, 6 and 8, the weak-coupling phase is different in 
three respects. First, the precisely linear dependence of (xx) seen in Figure || 
for iVt = 4, 6 and 8 as m^ai -^ becomes significantly non-linear for the 16^ 
lattice. Second, for P = 4.65 {xx) doubles as Nt is increased from 4 to 16. 

Finally for myaia = fngeaO' = 0.015 the hadron spectrum does not show 
the degree of parity doubling that might be expected for a phase in which 
chiral symmetry has been restored. For example, we see nearly exact parity 
doubling for larger P as is described in Section ^. In fact, the values of the 
hadron masses and (xx) seen here at /3 = 4.65 are more similar to those found 
for two flavors in the low temperature, chirally asymmetric phase using the 
same quark mass and lattice size. However, we do find parity doubling in the 
limit rriyai — * 0. We tentatively interpret these results as suggesting that these 
P = 4.65, 16'^ X 16 and x32 lattices lie in a transition region that occurs in 
the weaA;-coupling phase — a finite-temperature transition/cross-over leading 
to spontaneous chiral-symmetry violation on lattices of larger spatial and 
temporal extent. 



Let us now describe these results in greater detail. Tables [VIII| and |X 
give the particulars of the 16^ x 16 and 16^ x 32 calculations on which our 
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conclusions are based. Again, except where explicitly noted, we used the 
inexact R algorithm with the step size Ar = 0.0078125. Because of the 
limited length of some of these runs, the amount of data discarded as not in 
equilibrium will be discussed on a case-by-case basis below. 

Figure |] shows the evolution of (xx) starting from hot and cold starts 
for (3 = 4.65 on a 16^ x 16 lattice. The persistence of two phases over 
more than 300 time units, a time scale considerably greater than the initial 
thermalization time of < 50 time units, is evidence for a strong, first-order 



transition. As further evidence for a first-order transition. Figure |I0| shows 
an evolution for the 16^ x 32, f3 = 4.60 cold-start run in which a tunneling 
event occurs at r ;^ 250 time units, suggesting that the weak-coupling phase 
becomes unstable as (3 is decreased from 4.65 to 4.60. 

The critical coupling (3^ on the 16^ lattice is determined by the same 
procedure described earlier for A^^ = 4, 6 and 8. We created a "mixed" start 
by beginning with the weak-coupling configuration whose evolution is shown 
in the lower curve in Figure ^ and then varying j3 by hand for 55 time units 
to obtain a configuration with a value of (xx) lying midway between the 
strong- and weak-coupling values seen in the figure, i.e. {xx) ~ 0.2. This 
final configuration is then used as the beginning for the three different runs 



shown in Figure |TI]. This figure establishes (3c = 4.73(1), precisely the result 
found above for Nt = 8. Because of the significant time-step dependence seen 
earlier for (3c, we carried out this procedure a second time for Ar = 0.005 and 
determined for that case, (3c = 4.62(1), again in agreement with the Nt = 8 
result for that smaller time step. We conclude that this strong, first-order 
transition has become independent of the lattice size for Nf > 8 and hence 
is a T = transition. 

The lack of dependence of (3c on Nt seen for A^t > 8 is quite consistent 
with the At dependence of the discontinuity in the gauge action across the 
transition. For a normal, finite-temperature transition, an increase of A^^ 
by a factor of 2 would correspond to a decrease of the lattice spacing a by 
a factor of 2. Since for such a transition the discontinuity in the action is 
proportional to a physical latent heat, the jump in the action should decrease 
by a factor of 2^ = 16 when A^ increases from 8 to 16. In fact a similar factor 
of (4/6)'* is seen for the Nf = 4 latent heat when Nt is increased from 4 to 
6||16||. However, our 16^ x 8 and 16'* results given in Tables |V11| and^ show a 



large 20% jump in the action which changes relatively little between Nt 
and 16. In particular, for Nt = 8 the decrease in the action between (3 = A.70 
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and 4.75 is 0.1202(5) while for A^^ = 16 the difference between the action in 
the two metastable phases at /3 = 4.65 is 0.1150(2) — a decrease by 5% not 
by a factor of 16. 

Next let us consider the chiral-symnietry properties of the two phases sep- 
arated by this transition. We have computed {xx) ^^id the hadron spectrum 
in each phase for a number of valence quark masses. We work at /3 = 4.65, 
where there are two metastable phases as shown in Figure §. Although this 
choice of /? is below the critical value jSc = 4.73, at this P the weak-coupling 
phase shows no signs of instability either during the 342.5 time unit 16^ run 
(Figure ^ or the 865 time unit run on a 16^ x 32 lattice. These results, 
together with the /3 = 5.0 results described in the next section are given in 
Tables g, Q ^ and \Kin\ 



The second columns of Tables |KII| and pCIII| , for the strong-coupling phase, 
were obtained on a 16^ lattice from the microcanonical time range 210 — 380. 
The fitting was done by combining these results into blocks of 5 time units 
and the masses came from two parameter fits assuming a single propagating 
state. The tt mass was determined by fitting time separations 2 to 8 while 
the masses for the other states were obtained from time separations 1 to 5. 
The resulting vr masses appear to be well determined from the range of time 
separations available. However, the correlators for the other strong-coupling 
states decrease so rapidly with separation that useful information comes only 
from a few time separations. As a result, we are less certain that those masses 
have taken on truly asymptotic values. 



For the weak-coupling phase, given in the third column of Table pCll] and 



columns three through six of Table PCIII| we use the longer P = 4.65 run 
on a 16^ X 32 lattice discarding the first 382.5 time units for equilibration. 
With this larger time dimension, stable results are obtained for a larger 
number of masses. We follow a procedure to extract the masses similar to 



that used earlier [12]: the vr mass is obtained from a one-propagating-state, 
two parameter fit while the other masses come from a two-propagating-state, 
four parameter fit. We used a range of fitting separations, t^m < t < 16 as 
follows: For the quark masses 0.004, 0.01 and 0.015 and the A^ and A^' states 
we used tmm = 8 and for all other states tmm = 10. For the quark masses 
0.025 and 0.05 we used t^m = 10 for the vr state and tmin = 8 for the others. 
With these fitting ranges we obtain x^/dof < 2. All the masses given in 



Tables PCII| and PCIII] came from fits obtained by minimizing x computed 



from the full covariance matrix. Errors were determined with the jackknife 
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method using blocks of 15 time units (5 time units for the shorter 0.025 and 
0.05 runs) with corrections for autocorrelations in Monte Carlo time. 

As can be seen from Table |Xlij the strong- and weak-coupling phases 



have very different values for the hadron masses. The masses in the weak- 
coupling phase are lighter by about a factor of three, except the pion which 
we discuss in detail below. The vacuum expectation value {xx) ^-Iso becomes 
much smaller moving from strong to weak coupling. However, by itself, such 
a jump does not imply that the transition restores chiral symmetry. For 
example, naive scaling arguments would suggest that for small quark mass 
ixx) should decrease by a factor of 3^ when the hadron masses decrease by 
a factor of 3. 

To study the chiral symmetry of these two phases we show the linear 
extrapolations of {xx) ^"^^ "^^ as m^ai ^ in Figures |12] and |1^. The 
behavior of {xx) ^"^^ "^^ i^^ the strong-coupling phase (Figure |r^) is easily 
interpreted. We are seeing the usual consequences of spontaneous chiral- 
symmetry breaking — behavior quite similar to that shown for Nj = 2, f3 = 
5.7 in Figure |[|T2],g. 



However, the chiral properties of the weak-coupling, Nf = 8 phase are 



more ambiguous. Our earlier 16 x 32 results|12, 13| for Nf = 2, (3 = 5.7 
and TUseaO' = 0.015 are reproduced in column five of Table PCII| , allowing 
easy comparison of the chiral-symmetry breaking found in these two sets 
of spectra. Although the m^ — m„ and rup — uiai splittings are significant 
for the eight-flavor case, they are perhaps half the size of those seen in the 
Nf = 2 case. In the earlier Nf = 2 calculation we found {xx) = 0.0385(1) 
so the possible measure of chiral-symmetry breaking, {xx) l^pi is the same 
between the two calculations up to the 20% level. 

It is also of interest to ask how {xx) changes in the weak-coupling phase 
as lattice size increased. We can directly compare the Nt = A and Nt = 16 
results for (3 = 4.65 recognizing an increase from 0.0343(3) to 0.0711(4). 
Such an increase with increasing N^ might be interpreted as the onset of 
chiral-symmetry breaking as the temperature decreases for fixed jS. 

However, the notion that the weak-coupling phase seen for Nt = 16 and 32 
shows spontaneous symmetry breaking is not supported by the dependence 
of (xx) or the hadron spectrum on m^ai- As shown in Figure (xx) appears 
to approach as rriyai — ^ while the corresponding limit of 771^^ is small but 
non-zero. As can be seen, the small niyai -^ limit of (xx) is surprisingly non- 
linear for m.uaict as small as 0.004 but certainly appears to vanish. However, 
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the extrapolation of m^ is straight forward, 

mla^ = 0.0391(13) + 6.90(6)m„,rfa (x'^/dof = 1.9/3), (10) 

giving m'KO) ~ 30 standard deviations away from zero. This is in contrast 
with the behavior seen in the Nf = 2 weak-couphng phase shown in Figure ^ 
where, as is shown in Eq. (^, the extrapolated vr mass is consistent with 
zero. 



In Figure [T^ we show the masses of the three parity partners N' — N, 
Ai — p and cr — tt as a function of rriyai ■ Again even the considerable a — tt 
splitting disappears linearly as rriyai -^ 0. Quantitatively, fitting to the three 
smallest values of m^ai shown in the plot we obtain: 

m^a =0.218(5) + 10.7(4) m^^^a (xV^o/ = 2.9/1), 

m^a =0.220(4) + 16.5(4) m^„,a (xV^o/ = 5.0/1), 

ruptt =0.359(14) + 10.8(12) m^^^a {x'^ /dof = O.l/l), 
mA,a =0.371(30) + 14.1(24) m^„;a {x'^/dof = 0.003/1), 

m^a =0.704(35) + 11.2(26) m^^/a (xV^o/ = 0.01/1), 
mN,a =0.694(42) + 13.5(31) m^^/a (xV^o/ = 0.08/1), 



:ii) 



showing detailed chiral-symmetry restoration as m^ai -^ 0. Note the linear 
fit to m^ in Eq. (p!OD (using valence quark masses up to 0.050) and the linear 
fit to m^ in Eq. (pUD (using valence quark masses up to 0.015) give 0.198(3) 
and 0.218(5) respectively for m.„{0). 

We conclude that the strong-coupling phase seen at /? = 4.65 for Nt = 16 
and 32 shows clear spontaneous violation of chiral symmetry while the chiral 
symmetry of the weak-coupling phase is less obvious. However, the unusual 
rriyai dependence that we see in the weak-coupling phase for /5 = 4.65 is very 
much like the rrisea dependence found earlier in the four-flavor, 10^ x 6 work 
of DeTar and Kogut]^. Their results in the critical region [jS = 5.175) show 
a non-linear approach of (xx) to zero as rrisea —>' 0. Likewise, their tt and a 
screening lengths, while significantly non-degenerate for rrisea = 0.05, become 
equal when extrapolated linearly to rusea = 0. This behavior is precisely 
the niyai — > dependence that we see for these quantities. Therefore, we 
speculate that for Nt = 16 and (3 = 4.65, the weak-coupling phase is itself 
near a standard, finite-temperature transition region separating the chirally 
symmetric, weak-coupling behavior we see for Nt < 8 from a weak-coupling, 
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chirally asymmetric region that will be seen for Nt > 16 on significantly 
larger spatial volumes. This speculation is represented in Figure |l|, where the 
dashed line identifies a possible finite-temperature phase transition dividing 
the weak-coupling phase into low-temperature, chirally asymmetric (upper 
portion) and high-temperature, chirally symmetric (lower portion) phases. 
This line passes near Nt = 16 at (3 = 4.65 as is suggested by our results for 
these parameter values. 

5 High-temperature Region: Nf = 32 

In an attempt to understand the properties of the weak-coupling phase dis- 
cussed above for j3 ~ 4.65 on 16^ x 16 and x32 lattices, let us examine a 1325 
time unit calculation of (xx) ^"^^ hadron masses with (3 = 5.0 carried out 
on a 16^ X 32 lattice. As is discussed below, we find clear chirally symmet- 
ric behavior for this larger value of (3. The hadron screening lengths show 
complete parity doubling within errors, (xx) extrapolates linearly to zero as 
'fTT'vai -^ and rriTT varies little as rriyai -^ and has a relatively large rriyai = 
limit. 

These masses or screening lengths were determined from the evolution 
interval 605 — 1325. The results are shown in Tables LXllI and PCIVI. The 



fitting procedure is very similar to that used earlier: the vr-like states were 
determined from a two-parameter, single-state fit while the other states from 
a four-parameter, two-state fit. For all masses we used a fitting range from 
time separations 10 to 16. The quark masses of 0.01 and 0.015 were analyzed 
dividing the data into blocks of 15 time units while the shorter run with the 
valence mass of 0.004 used blocks of 5 time units. In contrast to the other 
mass fits discussed in this paper, the x^ values were very large, typically 
10 to 30 with 5 degrees of freedom. However, the jackknife errors for these 
X^ values were nearly as large as the x^ themselves and the x"^ computed 
ignoring off-diagonal terms in the correlation matrix are quite reasonable. 
We conclude that the fits are acceptable but that small poorly determined 
eigenvalues in the correlation matrix make the determination of x^ difficult. 
The hadron spectrum looks very much like that found in earlier calcula- 
tions in the plasma phaseQ. In particular the masses (or more accurately 
screening lengths) show remarkable parity doubling with the parity partners 
IT — a, p — Ai and N — N' having very nearly the same mass. Likewise we can 
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examine the extrapolation to zero valence mass of both rri^ and (xx) shown 



in Figure 15. Linear fits to the data yield 



mla^ =0.134(6) + 1.95(52) m„„;a {x^/dof = 2.0/1) 
{xx)a^ =0.00052(9) + 2.326(7) m^^/a {x^ / do f = 8.2/1) 



2/^.._«o/1^ (12) 



In marked contrast with the behavior seen for (3 = 4.65, mvr depends rather 
weakly on the valence quark mass, extrapolating to a value only 10% below 
the rriyai = 0.015 point, while (xx) extrapolates to a very small value. Given 
the non- vanishing of m^ as m^ai —>' 0, the statistically non-zero value of (xx) 
may reflect the use in Eq. (p^ ) of a fit neglecting the correlations between 
results for different valence quark masses. 

Although the behavior seen for jS = 5.0 is very clearly that expected from 
QCD at finite temperature, we should emphasize that our 16^ x 32 lattice 
is awkward to interpret as representing finite temperature. The nominal 
"temperature" direction with extent Nt = 32 and the required anti-periodic 
boundary conditions for the fermions is the longest dimension in the lattice. 
Probably the best interpretation of our space-time volume is as a 16^ x 32 
spatial volume with a temperature dimension corresponding to 16 lattice 
units. Clearly the behavior seen on this 16^ x 32 lattice may show significant 
finite-volume distortions relative to a proper, finite-temperature calculation 
on a Ng X 16 lattice with A^^ ^ 16. 

The contrast between the /? = 5.0 behavior just described and the P = 
4.65, weak-coupling phase discussed in Section ^ supports the hypothesis that 
for the 16^ lattice, /3 = 4.65 lies near a transition region. By increasing f3 from 
4.65 to 5.0 the degree of chiral symmetry has dramatically increased for fixed 
f^vai = 0.015 and the non-linear m^ai dependence of (xx) has disappeared. 

6 Possible A^j-Dependence of QCD 

In this section we present a possible picture of the Nf dependence of QCD 
that connects earlier work for Nf = 0, 2, 3 and 4 with the Nf = 8 results given 
here. Although the picture described below is supported by the presently 
available numerical results, it is far from unambiguously established by our 
current calculations. 

We would like to interpret the eight-flavor bulk transition seen here as 
an outgrowth of the strong- to weak-coupling cross-over region seen in pure 
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SU(3) gauge theory for (3 ~ 5.6. The variation seen in this region pro- 
vides a connection between strong couphng, where the scale of the physics is 
controlled by the lattice spacing, and weak coupling, where the scale is un- 
related to the lattice spacing. (The width of the cross-over region for SU(2) 
seen using the standard Wilson action can be altered by including an adjoint 
representation contribution to the action|T^.) As one passes through this 
region from strong to weak coupling, the hadronic energy scale (measured in 
lattice units) decreases at a rate faster than predicted by the perturbative 
renormalization group. This was seen quite clearly for T^ a in pure SU(3) by 
Kennedy et a/.[0. 

We hypothesize that adding additional light dynamical quarks to QCD 
promotes a rapid cross-over region to a phase transition, the first-order, bulk 
transition seen here, and that the addition of the quarks is similar to the 
effect of a non-zero adjoint action in the pure SU(2) case. The effect of this 
increasingly sharp cross-over region on the finite-temperature QCD phase 



transition might be deduced from Figure |Tg. Here we represent the cross- 
over region for a system of infinite space-time volume by the interval of (3 
between the vertical dotted lines. The more rapid variation of Nt with /3c 
within this region joins the relatively large value of Tc a for small (3 with a 
smaller value for large (3. If this region narrows as the number of flavors 
increases, sharpening into an actual discontinuity, the QCD phase transition 
for values of j3c within this cross-over region might be expected to sharpen as 
well. Such a sharpening of the transition as the number of flavors increases 
is certainly well established by current simulations^. 

Furthermore, such a narrowing of the cross-over region with increasing 
Nf would imply a corresponding increase in slope of A^j versus (3 within this 
region. In Figure ^ we plot the variation of Nt with (3c seen in simulations 
for zero|]T9|, two[^] and four||21| flavors together with that seen here for 



Nf = 8. For the two and four flavor cases, linear interpolation has been 
used to produce a value for Pc at niseaf^ = 0.015. The behavior predicted by 
the perturbative renormalization group is shown by the slopes of the dashed 
lines in the figure. This figure is consistent with the view that temperature 
dimensions between A''^ = 4 and 8 lie within this cross-over region and that 
the slope of the A^t- versus- /5c curve is increasing with Nf. 

Although far from well established, this picture is nicely consistent with 
the eight-flavor results presented in this paper. If the cross-over region shown 



in Figure |T6| shrinks to a vertical line as Nf -^ 8, the finite-temperature phase 
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transition effectively disappears for tlie corresponding interval of Nt, being 
engulfed there by the bulk transition. An eight-flavor phase diagram very 
much like that shown in Figure |I| results, in a manner that might be described 
as follows: 

1. For strong coupling and small values of Nt one expects to see a single 
phase transition that separates very different strong- and weak-coupling 
regimes. Nt should vary with (3 for small values of Nt characteristic of 
the small length scale important at strong-coupling (4 < A^^ < 8 in 
Figure ||) . 

2. For values of Nt larger than this strong-coupling length scale, the tran- 
sition becomes a bulk transition, with a fixed value of /3 = jSc- Now 
the dramatic change in hadronic length scale, which occurred rapidly 
in the cross-over region for Nf < 4, happens discontinuously across this 
bulk transition (a scale change by a factor of 3 in our case for A^t > 8). 
The finite-temperature, A't-dependent transition has disappeared. 

3. An apparently independent finite-temperature transition should occur 
in the weak-coupling phase at a much larger value of A^. This larger 
value of Nt (determined for /? near the bulk transition) should be related 
to the values of A^^ identified in 1 above. These two scales should be 
related by the same factor that describes the jump in the length scale 
of hadronic phenomena across the bulk transition. Thus, in our case 
we might expect a weak-coupling, finite-temperature phase transition 
to occur for j3 = 4.73 and values of Nt in the range of 3 (the jump 
in hadronic length scale) x8 (the Nt where the transition becomes A^ 
independent). In fact, as discussed in Section H, we have some evidence 
for such a weak-coupling, finite-temperature transition for /? = 4.65 and 
Nt ^ 16. A choice of /3 = 4.73 and A^ = 20 is used to locate the dashed 
curve in Figure |I|. 

We can compare the ratio Tc/rrip for this conjectured finite-temperature 
transition with the value for other numbers of flavors. For zero, two 
and four flavors, extrapolated to zero quark mass, Tc/rrip = 0.26, 0.19 
and 0.13||^. Using rup at /3 = 4.65 and assuming a monotonic decrease 
in Tc/mp with the number of flavors, we find Nt > 21. 
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Of course, a T = or bulk transition for which (3c becomes precisely 
independent of Nt should separate two phases each of whose properties are 
independent of A^^. This is not the case for a quark-gluon plasma, a natural 
candidate for the high temperature, weak-coupling phase represented by the 
lower right region of Figure 1. However, the 1/N^ behavior expected for the 
free energy of a quark-gluon plasma becomes sufficiently weak for A^j > 8 as 
to be completely consistent with the A^j-dependence of /3c that we see. In 
fact, /3c{Nt) for our four values of Nt is well fit by 

(3c{Nt) =4.737(6) - A0{3)/N^ xVdof = 0.6/2. (13) 

A final implication of our hypothesis, in analogy with the pure SU(2) 
case, is that for large Nt, the strong- and weak-coupling sides of the bulk 
transition may be continuously connected by using an action which includes 
single plaquette contributions from higher representations of SU(3). This is 
consistent with our picture that the bulk transition, for large enough lattices, 
is between two chirally asymmetric phases. 

7 Conclusion 

Our Nf = 8 studies are well summarized by the j3 — Nt phase diagram given 
in Figure |1]. Let us conclude with the following remarks: 

1. We have argued in Section | that the phase structure shown in Figure |I| 
may be quite consistent with the fiavor dependence of the QCD phase 
transition seen previously for Nf = 0, 2, 3 and 4. However, in that 
discussion we argued that the well-known strengthening of the finite- 
temperature transition for 4 < A'f < 8 that is seen with increasing Nf 
came from approaching the strong, Nf = 8 bulk transition. Thus this 
important feature, which dominates present lattice calculations, may be 
closely tied to a lattice artifact. The true, continuum, Nf dependence 
of the QCD phase transition may be quite different and may be seen 
only on much finer lattices. 

2. It is interesting to ask if such a bulk transition for eight fiavors may 
have already been anticipated. In fact there have been a number of 
papers that have explored possible phase structures for QCD with a 
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large number of quark flavors p3|, 0]. Surely interesting new behavior is 
to be expected as one approaches Nf = 16.5 where asymptotic freedom 
is lost. However, these theoretical studies typically predict a transition 
between a strong-coupling, chirally asymmetric phase with particle- 
like bound states and a zero-temperature, weak-coupling phase with 
Greens functions showing fractional anomalous dimension and lacking 
a particle interpretation. 

This behavior is not seen in the weak-coupling phase of Figure |1| for 
either f3 = 4.65 or 5.0. The masses or screening lengths given in Ta- 
bles pCli| , pClli| and pClV| are non-zero and come from fltting the correla- 



tion function to a function with exponential time dependence. Thus it 
is most natural to interpret these mass results as describing interacting, 
particle-like states with deflnite, non-zero energy eigenvalues. We be- 
lieve that the high-temperature, weak-coupling phase seen in our calcu- 
lations is quite conventional, very much like normal, high-temperature 
QCD seen for Nf < 4. 

3. A further argument for the flnite-temperature transition represented by 
the dashed line in Figure |I| is based on the 't Hooft anomaly conditions 1 22 1 



Since the 't Hooft anomaly conditions are inconsistent with a chirally 
symmetric eight-flavor, color SU{3) theory at zero temperature, we 
expect that the chirally symmetric phase to the right of the bulk tran- 
sition in Figure |l| cannot extend to zero temperature in the continuum 
limit. The flnite-temperature transition suggested by the dashed line 
in Figure ^ insures that this chirally symmetric phase is restricted to a 
region of non-zero temperature in the continuum, Nt —>■ oo limit. 

4. Our original objective in undertaking these eight-flavor calculations 
was to study the zero-temperature hadron spectrum for Nf = 8 to gain 
some quantitative insight into the effects of the fermion determinant 
in lattice QCD, hadron mass calculations. This objective has been 
frustrated by the existence of the Nf = 8 bulk transition. 

The decrease in Tc/rrip expected as Nf increases, requires that we in- 
crease nip a by lowering f3 or increase 1/Ta = Nf relative to calculations 
with smaller Nf. The Nf = 8 bulk transition prevents us from increas- 
ing rrip a by moving closer to the strong-coupling region and forces us to 
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work instead at larger Nt, in particular lattice sizes with spatial dimen- 
sions > A^i 3> 16. Thus it appears that both a proper demonstration of 
the finite-temperature, weak-coupling phase transition and such a low 
temperature study of hadron masses in eight-fiavor QCD forces the use 
of lattice volumes considerably larger than 16^ x 32. 

After the completion of this work, we became aware of a study of QCD 
with many fiavors of Wilson fermions in the strong coupling limit {f3 = 
0.0)|]23. Given the difference in the coupling and the type of fermions used. 



the overlap between [24] and the present work is unclear. 
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At 


Start 


P 


Total T 


Valence run parameters 


"^val 


Valence r 


0.0078125 


mixl 


4.50 


100 


0.004 


50-100 






4.55 


25 


- 


- 






4.58 


100 


0.004 


50-100 






4.60 


15 


- 


- 






4.65 


200 


0.004 


150-200 


niix2 


4.55 


50 


- 


- 






4.60 


50 


- 


- 






4.65 


25 


- 


- 



Table |I[ A list of the parameters for the runs with Nt = 4. The mixl start 
was produced by thermalizing a hot lattice at /5 = 4.5 for 80 time units 
and then varying f3 for 95 time units. The mix2 start was produced 
by thermalizing a cold lattice for 70 time units at /? = 4.5 and then 
varying f3 for 112.5 time units. 
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At 


Start 


P 


Total T 


Valence run parameters 


"^val 


Valence r 


0.0078125 


mix 


4.65 


25 


0.004 
0.010 


0-25 
0-25 


4.68 


50 


0.004 
0.010 


0-50 
0-50 


4.70 


50 


0.004 
0.010 


0-50 
0-50 


4.73 


50 


0.004 
0.010 


0-50 
0-50 


cold 


4.70 


100 


0.004 
0.010 


0-100 
0-100 


hot 


4.70 


100 


0.004 
0.010 


0-100 
0-100 



Table J^. A list of the parameters for the runs with A^^ = 6. The mixed 
start was produced by thermalizing a hot lattice for 50 time units at 
/3 = 4.70 and then varying (3 for 135 time units. 
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At 


Start 


/3 


Total T 


Valence run parameters 


"^val 


Valence r 


0.0125 


mix 


5.20 


25 


0.004 


0-25 






5.25 


50 


0.004 


0-50 






5.28 


50 


0.004 


0-50 






5.30 


50 


0.004 


0-50 






5.35 


50 


0.004 


0-50 


0.0078125 


mix 


4.70 


100 


0.004 


25-100 






4.73 


125 


0.004 


75-125 






4.75 


165 


0.004 


40-165 






4.80 


50 


0.004 


25-50 


0.005 


mix 


4.60 


93.75 


0.004 


0-93.75 






4.63 


62.5 


0.004 


0-62.5 






4.65 


62.5 


0.004 


0-62.5 






4.73 


93.75 


0.004 


0-93.75 


0.002 


cold 


4.59 


50 


0.004 


0-50 


mix 


4.55 


105 


0.004 


0-105 






4.58 


25 


0.004 


0-25 






4.60 


25 


0.004 


0-25 






4.62 


145 


0.004 


0-145 


exact 


mix 


4.58 


75 


0.004 


0-75 






4.60 


75 


0.004 


0-75 



Table [II 



A list of the parameters for the runs with Nt = 8. The mixed 
start was produced by thermalizing a cold lattice for 40 time units at 
/3 = 4.60 with At = 0.0078125 and then varying (3 for 60 time units. 
The At = 0.005 runs had a trajectory length of 0.625 time units. 
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Nt 


Ar 


Pc 


4 


0.0078125 


4.58(1) 


6 


0.0078125 


4.71(1) 


8 


0.0125 
0.0078125 
0.005 
0.002 

exact 


5.29(1) 
4.73(1) 
4.64(1) 
4.59(1) 
4.59(1) 


16 


0.0078125 
0.005 


4.73(1) 
4.62(1) 



Table [V. Values for (3^ for Nt = 4, 6, 8 and 16. 



Ar 


Start 


P 


''"eq 


Action 


ixx) 


Valence results 


"^val 


(XX)val 


0.0078125 


mixl 


4.50 
4.65 


50 
50 


0.6344(6) 
0.4939(1) 


0.418(2) 
0.0343(3) 


0.004 
0.004 


0.410(4) 
0.0094(4) 



Table |V]. Results for (xx) ^^^ the gauge action for Nt = 4. 



Ar 


Start 


P 


''"eq 


Action 


ixx) 


Valence results 


"^val 


(XX)val 


0.0078125 


hot 


4.70 


50 


0.5995(4) 


0.378(2) 


0.004 
0.010 


0.371(2) 
0.375(2) 


cold 


4.70 


50 


0.4867(2) 


0.0445(4) 


0.004 
0.010 


0.0126(4) 
0.0302(4) 



Table [V^. Results for (xx) and the gauge action for A^i = 6. 
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At 


Start 


^ 


''"eq 


Action 


{xx) 


Valence results 


"^val 


(XX)val 


0.0125 


mix 


5.25 
5.35 


30 
30 


0.6908(3) 
0.4168(4) 


0.472(3) 
0.0272(5) 


0.004 
0.004 


0.467(5) 
0.0077(5) 


0.0078125 


mix 


4.70 
4.75 


70 
100 


0.6006(3) 
0.4804(4) 


0.378(1) 
0.0465(3) 


0.004 
0.004 


0.372(2) 
0.0131(3) 


0.005 


mix 


4.60 
4.73 


62.5 
62.5 


0.5793(5) 
0.4838(1) 


0.348(2) 
0.0493(6) 


0.004 
0.004 


0.341(3) 
0.0132(6) 



Table [VII| . Results for {xx) and the gauge action for Nt = 8. 
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At 


Start 


P 


Total T 


Valence run parameters 


"^val 


Valence r 


0.0078125 


cold 


4.65 


342.5 


0.004 
0.010 


0-215 
215-342.5 


hot 


4.65 


380 


0.004 
0.010 


0-210 
210-380 


5.00 


40 


0.010 


0-40 


mix 


4.71 
4.73 
4.75 


50 
100 
100 


0.004 
0.004 


50-100 
50-100 


0.005 


cold 


4.62 


350 


0.004 


0-350 


hot 


4.62 


100 


0.004 


0-100 


mix 


4.60 
4.62 
4.63 
4.64 


50 
25 
25 
50 


0.004 
0.004 
0.004 
0.004 


0-50 
0-25 
0-25 
0-50 



Table |VllJj . A list of the parameters for the runs with A^^ = 16. The mix 
start was produced by thermalizing a cold lattice for 342.5 time units 
at /3 = 4.65 with At = 0.0078125 and then varying (3 for 55 time units. 
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At 


Start 


P 


Total T 


Valence run parameters 


"^val 


Valence r 


0.0078125 


cold 


4.60 


485 


0.004 
0.010 


0-355 

355-485 


4.65 


865 


0.004 


232.5-382.5 
532.5-765 


0.010 


0-232.5 

382.5-765 


0.025 
0.050 
0.100 
0.200 


532.5-865 

765-865 

765-817.5 

817.5-865 


4.70 


100 


0.004 


0-100 


4.80 


75 


0.004 


0-75 


4.90 


100 


0.004 


0-100 


5.00 


1325 


0.004 
0.010 


1182.5-1325 
397.5-1182.5 



Table [X 



A list of the parameters for the runs with A^^ 
run tunneled around r = 250. 



32. The p = 4.60 



Ar 


Start 


P 


''"cq 


Action 


ixx) 


Valence results 


"^val 


(XX)val 


0.0078125 


cold 


4.65 


100 


0.4961(1) 


0.0711(4) 


0.004 
0.010 


0.0259(5) 
0.0541(5) 


hot 


4.65 


100 


0.6111(2) 


0.3953(4) 


0.004 
0.010 


0.3900(11) 
0.3933(6) 


0.005 


cold 


4.62 


100 


0.5006(1) 


0.0829(6) 


0.004 


0.0399(8) 


hot 


4.62 


75 


0.5631(3) 


0.314(1) 


0.004 


0.305(2) 



Table ^. Results for (xx) ^^id the action for Aj = 16. 
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At 


Start 


(3 


Tcq 


Action 


ixx) 


Valence results 


"^val 


(XX)val 


0.0078125 


cold 


4.60 


325 


0.62031(6) 


0.4075(3) 


0.004 
0.010 


0.4038(16) 
0.4053(4) 


4.65 


382.5 


0.49552(6) 


0.0687(1) 


0.004 
0.010 
0.025 
0.050 
0.100 
0.200 


0.0234(3) 
0.0506(3) 
0.0983(1) 
0.1532(2) 
0.2255(2) 
0.3390(2) 


5.00 


250 


0.45125(2) 


0.03539(3) 


0.004 
0.010 


0.00965(9) 
0.02386(4) 



Table ^T. Results for {xx) ^-^d the action for Nt 
tunneled at about r = 250. 



32. The p = 4.60 run 



iJ"') 


/3 = 4.65 strong. 


/? = 4.65 weak 


P = 5.00 


/? = 5.7 


7r(0-+) 


0.297(1) 


0.378(2) 


0.405(3) 


0.293(2) 


^2(0-+) 


1.60(8) 


0.471(7) 


0.434(10) 


0.333(3) 


a(0++) 


- 


0.465(3) 


0.415(4) 


0.487(12) 


P(l-) 


1.41(1) 


0.522(7) 


0.484(7) 


0.455(8) 


P2{1-') 


1.64(2) 


0.521(4) 


0.490(5) 


0.452(7) 


^l(l+-) 


- 


0.582(11) 


0.491(7) 


0.594(22) 


N{f) 


2.29(7) 


0.872(10) 


0.807(7) 


0.685(10) 


N'Cf) 


- 


0.896(13) 


0.810(6) 


0.833(38) 


5i(l++) 


- 


0.586(26) 


0.512(8) 


0.596(28) 



Table Um 



Hadron masses for quark mass m-uaid = rngeaO' = 0.015. The 
second column was obtained on a 16"^ lattice beginning with a hot 
start and the third and fourth on 16^ x 32 lattices with a cold start. 
For reference, the right column lists the results of our earlier Nf = 2, 
16^ X 32 calculation with m.^n' 



TUyaia = 0.015 11, 12, 13 
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(J^-') 


/5 = 4.65 strong 


/5 = 4.65 weak 


7r(0-+) 


"^m/O = 0.01 


= 0.004 


= 0.01 


= 0.025 


= 0.05 


0.243(1) 


0.257(4) 


0.328(2) 


0.465(4) 


0.619(2) 


7r2(0-+) 


1.65(12) 


0.369(32) 


0.420(6) 


0.564(5) 


0.767(7) 


a(0++) 


- 


0.284(3) 


0.391(3) 


0.587(5) 


0.790(6) 


P(l--) 


1.40(2) 


0.404(11) 


0.465(9) 


0.632(6) 


0.813(4) 


P2(l--) 


1.63(2) 


0.407(11) 


0.468(7) 


0.634(7) 


0.817(7) 


^l(l+-) 


- 


0.428(28) 


0.511(14) 


0.751(40) 


0.951(26) 


N{f) 


2.40(15) 


0.747(30) 


0.818(21) 


0.992(15) 


1.246(8) 


N'Cf) 


- 


0.744(34) 


0.836(28) 


1.043(31) 


1.304(29) 


i?i(l++) 


- 


0.429(21) 


0.511(24) 


0.765(45) 


1.045(73) 



Table pCllij . Hadron masses calculated with a variety of valence quark 
masses. The masses quoted in the second column were obtained on 
a 16^ lattice beginning with a hot start while those in columns three 
through six came from a cold start using a 16^ x 32 lattice. 



(J^-') 


p= 


5.00 


7r(0-+) 


m^aiO' = 0.004 


rriyaia = 0.01 


0.386(9) 


0.389(3) 


7r2(0-+) 


0.500(108) 


0.422(27) 


a(0++) 


0.384(9) 


0.394(5) 


P(l--) 


0.474(6) 


0.465(10) 


P2(l-) 


0.485(8) 


0.481(8) 


^i(l+-) 


0.475(6) 


0.469(10) 


iv(r) 


0.783(22) 


0.785(8) 


N'Cf) 


0.789(28) 


0.786(7) 


5i(i++) 


0.487(8) 


0.491(11) 



Table PC1V| . Hadron masses for two different valence quark masses on a 
16^ X 32 lattice. 
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Figure Captions 



Figure 1: A phase diagram in the j3 — Nt plane for eight-flavor QCD in 
infinite spatial volume consistent with the results presented here. Nt is the 
temporal extent of the lattice and (3 = 6/g'^ is the lattice-coupling strength. 
The solid line, becoming vertical for Nt > 8 locates a "zero-temperature", 
first-order transition — a lattice artifact. The dashed hne suggests a possible, 
continuum finite-temperature transition that occurs in the weak-coupling 
phase. The system shows chiral symmetry to the right of and below this 
dashed line while we speculate that chiral symmetry will be spontaneously 
broken to the left of and above this line. The solid squares label parameter 
values where we have performed simulations, while the open squares locate 
critical values. 



Figure 2: The quark mass dependence seen for (xx) ^^id m^ in earlier two- 
flavor calculations. The ma = 0.01, 0.015, 0.02 and 0.025 points are cal- 
culations properly including the effects of dynamical quarks p^ , while the 
ma = 0.004 point is obtained using that value in the explicit quark prop- 
agators but the value ma = 0.01 in the fermion determinant. The lines 
correspond to the fits in Eq. (||). 



Figure 3: Generation of a mixed-phase configuration from a hot start for 
Nt = 4. The lower, cold-start trajectory establishes the value of (xx) fo^^ 
the chirally symmetric phase while the upper trajectory both gives a value 
of (xx) ill the symmetry broken phase and with a subsequent tuning of /3 
becomes our candidate "mixed" phase. 
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Figure 4: Determining (3c for a mixed-phase configuration generated from 
a liot starting lattice witli Nt = 4. From tliis figure we conclude that the 
critical value of f3 is (3c = 4.58(1). 

Figure 5: Determining (3c for a mixed-phase configuration generated from 
a cold starting lattice with Nt = 4. From this figure we conclude that the 
critical value of P is Pc = 4.58(1). 

Figure 6: Determining (3c for a mixed-phase configuration generated from 
a cold starting lattice with Nt = 8. From this figure we conclude that the 
critical value of (3 is (3c = 4.73(1). 

Figure 7: (3c versus Ar for Nt = 8. The curve is a quadratic fit to the points 
At = 0.002, 0.005 and 0.0078125. 

Figure 8: Extrapolation of (xx) for A^t = 4, 6 and 8 as a function of valence 
quark mass, m^ai- The fits are forced through the origin. 

Figure 9: The evolution of {xx) fo^^ two independent Monte Carlo runs on a 
16^ lattice at (3 = 4.65. The upper curve represents a run begun with a hot 
start while the lower curve began with a cold start. 

Figure 10: The evolution of (xx) from a cold start on a 16^ x 32 lattice with 
(3 = 4.60. We interpret the jump seen at r ^ 250 as tunneling from the 
metastable, weak-coupling phase to the stable, strong-coupling phase. 

Figure 11: Determining (3c for a mixed-phase configuration generated from 
a cold starting lattice with Nt = 16. From this figure we conclude that the 
critical value of /5 is Pc = 4.73(1). 

Figure 12: Linear fits to (xx) ^-^id m'^ in the strong-coupling phase at /3 = 
4.65 on a 16^ lattice. The m'^ fit is forced through the origin and has x^ /dof = 
0.1/1. 
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Figure 13: (xx) ^^^d m^ in the weak-coupling phase at /3 = 4.65 on a 16^ x 32 
lattice. The line is a linear fit to m^. 



Figure 14: Values m^, rrifj, nip, tjiai, ttin and mAr/ plotted versus m^aiO' for a 
16^ X 32 lattice in the weak-coupling phase with /? = 4.65. The lines shown 
correspond to the fits in Eq. ([TT|). 



Figure 15: (xx) ^^'^ '^^ plotted as a function of m^aiO' for /3 = 5.00 on 
a 16'^ X 32 lattice. The straight lines shown are least squares fits to the 
three mass values m^aiO' = 0.004, 0.01 and 0.015. Both quantities show the 
behavior expected in a chirally symmetric phase. 



Figure 16: A sketch of the jSc dependence of Nt for the finite-temperature 
QCD phase transition for Nf < 4. The area between the vertical, dotted 
lines represents the region of (3 for which infinite space-time volume systems 
show cross-over behavior from strong to weak coupling. The increased slope 
of the A'"t-versus-/3c curve in this region displays the well-established non- 
scaling behavior seen for Tea in this region. The dashed line has a slope 
predicted by the perturbative renormalization group. 



Figure 17: Values of Nt versus /3c are plotted for zero|]T9l, two|l20[|, four [pi] 



and eight flavors. The dashed lines have slopes predicted by the perturbative 
renormalization group and have been located to show the possible weak- 
coupling behavior of the adjacent curve. 



36 



